How Colors Influence Numbers: Photon Statistics of Parametric Downconversion 
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Parametric downconversion (PDC) is a technique of ubiquitous experimental significance in the 
production of non-classical, photon-number correlated twin beams. Standard theory of PDC as a 
two-mode squeezing process predicts and homodyne measurements observe a thermal photon number 
distribution per beam. Recent experiments have obtained conflicting distributions. In this paper, we 
explain the observation by an a-priori theoretical model solely based on directly accessible physical 
quantities. We compare our predictions with experimental data and find excellent agreement. 

PACS numbers: 42.50.Ar 89.60.Gg 



Introduction Spectral properties of states gener- 
ated by x*-^-* nonlinearities are traditionally studied using 
homodyne detection. Unfortunately, this standard tech- 
nique implicitly restricts the observation to an effective 
single spectral mode imposed by the single local oscil- 
lator. Avalanche photo diodes (APDs) [1], in contrast, 
are sensitive on all modes generated by sources of cur- 
rent experimental significance, and uncover richer spec- 
tral properties. This sub-structure is currently usually 
neglected or only treated effectively, although it impacts 
security proofs of quantum key distribution or the valid- 
ity of fundamental quantum measurements, for example. 

In this paper, we present an a-priori theoretical expla- 
nation that connects the spectral structure of PDC states 
with the photon number distribution (PND), which is a 
commonly employed resource. Recent experiments have 
observed that the PND for multi-mode sources differs 
markedly from the prediction of the single-mode standard 
model Our approach explains this behavior by de- 
composingthe state into a set of independent two-mode 
squeezers [1, Q akin, but not completely identical to the 
Bloch- Messiah decomposition. The PND is inferred from 
the well-known properties of these independent contribu- 
tions. In contrast to previous efforts [![ [1| , our approach 
is the first to enable, to our knowledge, the quantita- 
tive computation of photon number statistics without as- 
sumptions or fitting of non-physical parameters. This is 
important for a wide class of experiments ranging from 
fundamental to highly applied because they require a 
complete understanding of the internal structure of PDC 
states to fully exploit their quantum features. 



Decomposition A multi-mode type-II down- 
conversion process is most conveniently stud- 
ied using the interaction Hamiltonian iJint(i) = 



H.c. 0, 

where the subscripts denote pump, signal, and idler, 
respectively, and the tensor x*^^^ represents the second- 
order nonlinear susceptibility. By assuming a classical 
pump and a frequency-independent x^^^ ™- the spec- 
tral range of interest, it can be shown [s'] that with 



jy d?xx^^)Ei+^ (f , t)E^s-^ {x, t)E[-\x, t) 



Hi = J^^dt'H,,,tit'), 

Hi = C 11 dcc;idcj2/(wi,c^2)a^K)6^(t^2) +H.C., (1) 



where a){LOi) and 6^ (0^2) are field operators that create 
a monochromatic photon with frequency Ui in the signal 
and idler modes a and b. f{L0i,uj2) is the spectral distri- 
bution function (SDF) of the single photon contribution, 
and C — C {x'''^\ \/lp) is a coupling constant that de- 
pends on the strength x^^'' of the nonlinear susceptibility 
and on the pump intensity 0, 0| . The time-propagated 
state is computed by \ip) = Texp{{ih)~^Hi)\^{to)), 
where we assume that the pulse has completely left the 
crystal and the interaction is finished. Following [sj, the 
time-ordering T can be omitted because the Hamiltonian 
approximately commutes with itself at different times 
and the corrections are therefore negligible. 

To express Hj in a more convenient form, we use the 
Schmidt decomposition, uniquely defined by 



JV-l 



(2) 



ri=0 



where the Schmidt modes {^n\'-^i)} and {^i^' (^2)} are 
two sets of orthonormal bases with respect to the in- 
ner product, and the Schmidt eigenvalues A„ are real ex- 
pansion coefficients that satisfy A„ = 1. The salient 
feature of Eq. ([2]) is that only a single summation index 
is required, and not two as for a regular change of basis. 
The decomposition is guaranteed to exist for a large class 
of systems under very general assumptions ^1^ . For sim- 
ple systems that require only a few Schmidt modes (i.e., 
N is small), the decomposition can be numerically com- 
puted by solving a set of coupled integral equations [ll[ . 
For systems that require a large N, it is usually easier 
to perform a singular value decomposition (SVD), see 
Ref. and below for more details. 

We define effective single-mode field operators (some- 
times also called pseudo-boson operators) by 



(3) 



2 



and similarly for i3„. Because the spectral distribution 
functions are orthonormal, that is, = Sij, it is easy 

to verify that the operators fulfill the canonical commu- 
tation relations = iSjk and [Aj, A^] = 0. More 
details about this notation are provided by Ref. fl^ . 

By rewriting Hj in Eq. ([1]) using the Schmidt decompo- 
sition ([2]) for f{uji,uj2) and the definition of pseudo-boson 
operators in Eq. ([3]), we obtain 



1 



n—0 



(4) 



The two-mode squeezing operator for spectral effec- 
tive single modes A, B is defined by 5'AB('?n) = 
exp(-77„iti?t + 77;iS), where 77„ = CVX^/iifi) = 
r„expz(/5„ is a complex number. Because = 

for j ^ k, the state after the interaction is a tensor prod- 
uct of independent two-mode squeezers [2lj: 



N-l 



SA„B,SVn)\lp{to))- 



(5) 



Notice that it follows from this decomposition that the 
SDF is identical for all orders of photon number contribu- 
tions because creation operators that belong to different 
distribution functions are never mixed. 22 1 



Computing Statistical Distributions For two- 
mode squeezed states, the PND in each mode is thermal, 
that is, for the state 



w = SAB{v)m 



[ Kn\n,n), 



(6) 



the distribution is given by p{n) — |k„P — 
sech^ r tanh^" r for one output mode, that is, = 1. 
Consequently, the photon number distribution of the 
multi-mode state ([5|) is given by the convolution of the 
distributions of all independent squeezers. Assume that 
p^j. (n) denotes the PND of the fcth squeezer with spectral 

modes ^1 . The overall PND is then given by 



N-l 



(7) 



&£nhN m=0 



where n\- N denotes the set of all partitions of n into N 
parts. The distribution p^{n) is consequently the convo- 
lution of all probability distributions p^i{n). 

Two special cases follow directly from Eq. ([7]): When 
only a single effective mode contributes {N = 1), the 
resulting distribution exhibits thermal behavior. When 
the physical process requires a very large number of effec- 
tive modes (N oo), the resulting PND is Poissonian, 
because it is known that a convolution of thermal dis- 
tributions converges to a Poissonian distribution in this 
limit [141] . 



Computing the convolution in Eq. ([7]) involves sum- 
ming over numerous contributions. This is consider- 
ably simplified by using generating functions. For coeffi- 
cients they are given by the formal power series 14 1 
g{C) = SnP('^)C"- The individual coefficients can be re- 
covered via p{n) = ^-^pr 5(01^=0- ^'^^ thermal dis- 
tribution of a two-mode squeezer, the series converges an- 

1 2 

alytically to gk{C) = i-<;°tanhV,- ' where rt is the strength 
of the fcth squeezer. The generating function for a convo- 
lution of N thermal distributions is J^^^^ fffclO? the 
resulting photon number distribution is consequently 



1 / a" ' 



(8) 



Let us now turn our attention to an example illustrat- 
ing our considerations. Assume that the SDF is given 
by a two-dimensional, real-valued Gaussian distribution 
(this is not a restriction because the methods also work 
for complex, non-Gaussian SDFs). This approximation is 
commonly used 11, to provide a convenient param- 
eterization of type-II PDC processes. Especially, it is 
possible to perform an analytical Schmidt decomposition 
(a similar approach is used, for instance, in Ref. .15]). 
We use the parameters ct^ and ay to specify the spectral 
widths of signal and idler, while 9 denotes the rotation 
with respect to the x axis. This form is illustrated in 
Figure [31 

Let us choose = 25 and cjy — 1, which are the pa- 
rameters depicted in the inset of Figure [3l The Schmidt 
number K = 1/ 'Y^^ \\ is computed from the eigenval- 
ues A„ of the Schmidt decomposition. It is a measure for 
the number of effectively contributing spectral modes and 
thus of inherent spectral correlations of the physical pro- 
cess [ll| (notice that we could have also considered an 
entanglement monotone like the logarithmic negativity 
for this purpose). For 6* = 0, the state exhibits no spec- 
tral correlations, and a single Schmidt mode suffices for 
the decomposition. By rotating the SDF from 6* = to 
9 = 7r/2, the correlations increase to their maximal value 
at 9 = 7r/4, and decrease again until the SDF becomes 
separable for 9 — Tr/2. This implies thermal statistics for 
9 = and 9 = tt/2, and maximal similarity to Poissonian 
statistics for 9 = 7r/4. The coupling and pump intensity 
are, for better comparability, chosen such that n = 1 for 
all PNDs. Figure [1] illustrates the arising distributions. 

To quantify the difference between convoluted and 
Poissonian or thermal distributions, we employ the vari- 
ational distance defined for two probability distributions 
pi, p2 as Apj^p2 = J2n -P2in)\. Two distributions 

are completely identical if and only if A = 0. Figure [2] 
compares the difference of the convoluted distribution to 
the above-mentioned special cases for a growing Schmidt 
number K, that is, a growing number of Schmidt modes 
achieved by rotating the Gaussian SDF for 6' = to 
9 = 7r/4. 
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FIG. 1: Photon number distribution depending on the num- 
ber of effectively contributing modes as given by the Schmidt 
number K (and thus on the angle of the SDF) of a type-II 
PDC process. The x axes depict photon numbers, whereas 
the y axes show probabilities. 




FIG. 2: Solid line and dashed line show the distance between 
the convoluted photon-number distribution and Poissonian or 
thermal statistics, respectively, plotted against Schmidt num- 
ber. For a single effective mode, the distribution is exactly 
thermal, but the more modes contribute, the closer it gets to 
a Poissonian distribution. The inset fixes 6 — tt/A and varies 
CT^, which is drawn on a logarithmic scale. 



Once again, we emphasize that the shift towards a 
Poissonian distribution is inherent in the physical pro- 
cess and not caused by any experimental imperfections. 

Comparison with Experimental Data We have 
also performed a comparison of experimentally measured 
photon number statistics with the predictions of our the- 
ory. A photon- number resolving fiber-loop detector [ij in 
combination with highly efficient waveguides was used to 
record the distribution. The detection method is resilient 
against loss and allows us to eliminate the correspond- 
ing effects when ensemble measurements are performed. 



Ref. [I4I shows the experimental details of state genera- 
tion, and describes the measurement procedure. Fig- 
ure [3] compares the experimentally observed distribution 
with the theoretical prediction at various pump powers. 
As is immediately obvious from the figure, they are in 
excellent agreement. 



Photon Number Distribution 




Ptioton number 

FIG. 3: (Color online) Comparison between experimentally 
measured and theoretically obtained photon number distri- 
butions for a multi-mode PDC process at various pump 
strengths. The bottom inset shows the real part of the joint 
spectral intensity, while the top inset demonstrates the pa- 
rameterization of the analytical Gaussian approximation of 
the SDF. Loss inversion and error estimation was performed 
using non-negative least squares optimization. 

To avoid the necessity of fitting any effective param- 
eters, we have obtained an exact numerical decomposi- 
tion using SVD techniques. After discretizing the SDF 
on a grid Mmn of size 1500 x 1500, the matrix is de- 
composed as M = UT,V\ where U, V are unitaries and 
E = diag(\/A7, . . . , \/X/v) is a real diagonal matrix [3]. 
Extensive checks that the decomposition converges (and 
also con verg es to the proper value) have been performed. 



see Ref. [17| for details. 



Notice that the decomposition of the spectral distribu- 
tion does not depend on the pump intensity, which means 
that the composition {A„} of the PND is fixed for the 
physical process. However, the observed mean value of 
the PND does depend on the pump intensity, and Fig. [3] 
shows a shift toward larger mean photon numbers for 
larger pump intensities as expected. 

For higher pump powers, photon-number resolved de- 
tection is not possible anymore. To check the theory in 
this regime, we have used a set of mean photon number 
(h) measurements instead. The coupling constant C as 
defined in Eq. [1] can be inferred from the decomposed 
SDF for each n for a given pump power by a numerical 
optimization process[23j . The result is shown in FigurelH 
Again, very good agreement between theory and experi- 
ment is achieved. 
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FIG. 4: Relation between pump power and coupling param- 
eter. The expected square-root dependency [Oj is correctly 
obtained for a wider range of pump intensities than can be 
resolved with current TMDs, which ensures the validity of our 
approach also for high powers. Notice that this knowledge 
would also allow for computing the expected mean photon 
number for a given pump power or a determination of x^'^\ 
as described in [laj . 



Conclusions We have shown how to decompose a 
multi-mode PDC process into independent two-mode 
squeezers operating on effective single modes, and how 
this explains why the photon number distribution of the 
process can exhibit any form ranging from purely thermal 
to purely Poissonian. We have underlined the validity of 
the theory by comparing the predictions to an experi- 
mentally measured photon number distribution. Addi- 
tionally, we have compared theory and experiment for 
larger pump powers. 

Appendix A two-dimensional Gaussian distribution 
in a suitable parameterization is given by 



1 



: exp(— aa;^ — 2hxy — cy^), 



aie, a^,c7y) = cos^ e/{2al) + sin^ 0/{2al), 

-sin26'/(4f72) +sin26'/(4a^), 



'V 



sin^ 0/{2al)+coa^ 0/{2al) 



Without getting into details of the alge- 
bra involved, we remark that by starting from 

Mehler's formula ^ T.n=o Hn{x)Hn{y)^-^ = 

_L=exp(-2:!£!^fl^¥±2Vi (i7„(x) denotes the 

Hermite polynomial of nth order), it is possible to bring 



Xnfi'\x)fi,^>{y) 
(i) 



(2), 



1+7 



2n 



ac Uj^cTy 



f{x,y) into the form f{x,y) = Y,°^= 
The coefficients A„ are given by A„ 
where 7 — ~2v^+V4ac-4fc'' ^ is^Yace the set {A„} contains 
all information required for our calculations, the exact 
form of^*''(-) is not of interest here, but can be found 
in Ref. [H. 
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SAjBj iv) = Scj (?;)<8>>S'dj (»?), that is, a product of two in- 
dependent effective single-mode squeezers [l^. Using this 
transformation, we obtain the standard Bloch-Messiah 
decomposition (see, e.g., Refs. 0,0]). Owing to the cou- 
pling of the signal mode Aj with the idler mode Bj , this 
form delivers the joint photon number distribution for 
signal and idler, pjoint (n) . It is connected to our distribu- 
tion via Pjoint (2n) = p{n), and pjoint(2n + 1) = 0. When a 
degenerate PDC process (including type-I) with Aj — Bj 
is considered, the decomposition in Eq. (j4| automatically 
leads to the Bloch-Messiah decomposition. 
This property can already be inferred from Eq.|4]by defin- 
ing the operators Xii^'' = Aj^Bj^, K^~^ = A„B„, and 
Kn'' = l/2(yiJjA„ -I- Bl^Bn) that share the commutation 



relations [Kli''' . 



±i:i*' and [K'„ 



(-) 



of an su(l, 1) Lie algebra, which allows us to apply a 
specific exponential operator disentangling formula [l9l . 
A5.18] from which the desired property is readily derived. 
If x'^' were known precisely, the optimization would not 
be necessary 



